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Abstract 


The presented paper shows a possibility of creation of a useful and efficient model of conductivity for composite polymeric electrolytes. It deals 
with a more elaborate and versatile mesoscale approach (in comparison to the already applied Effective Medium Theory approach for which each 
new composite microstructure forces a new mixing rule development) based on the Random Resistor Network method. The introduced model can 
predict the ac conductivity values for a whole spectrum of frequencies for different composites with a wide range of possible microstructures. It is 
based on complex current iteration procedure for a randomly generated three-dimensional network of impedances constituting a digital image of 
the electric properties of the virtual sample studied. A spatial unit of the material present in the sample is defined as a parallel RC connection. 
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1. Introduction 


Semicrystalline polymeric electrolytes based on poly- 
(ethylene oxide) matrix reveal good mechanical and electro- 
chemical properties unfortunately together with low electrical 
conductivity being the main limitation of their application. An 
improvement of electrical conductivity with a parallel further 
increase of the stability can be achieved by the formation of 
composite materials. They are produced through the addition 
of an electrically insulating filler (ceramic powders, immiscible 
polymers) to the pristine solid electrolyte. The obtained mate- 
rials can be applied to various electrochemical systems such 
as lithium—polymer batteries [1]. The micromorphology of the 
sample can be described as a suspension of the randomly dis- 
tributed filler grains in the polymer matrix. An explanation of 
the conductivity enhancement phenomenon [2,3] which is now 
widely accepted for the solid samples is based on the assumption 
that the formation of amorphous domains on the matrix-filler 
interface occurs due to the polymer crystallization process inhi- 
bition on the grain surface. As the amorphous phase conductivity 
can be assumed on the basis of the NMR data to be about three 
orders of magnitude higher than the crystalline one, continuous 
percolation paths can be formed in the sample. A wide range of 
examples [3—6] can be found in literature to confirm the fact that 
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the addition of nonconductive ceramic fillers such as alumina or 
silica to an electrolyte based on the PEO matrix leads to the con- 
ductivity improvement. A peak value in ambient temperature is 
near to 1075 Scm7! and can be obtained for composites con- 
taining 10-15% by weight of the dispersed particles. This value 
is about three orders of magnitude higher than for the pristine 
electrolyte and reaches the value characteristic for the amor- 
phous phase of the PEO. For higher grain concentrations one can 
observe a decrease of conductivity as an effect of the synergetic 
combination of the polymer material dilution and amorphous 
phase stiffening. A semi-quantitative dependence of conductiv- 
ity of the composite on the amount of the filler added can be 
divided into three separated regions. In the small filler concen- 
tration range only a very slight conductivity increase can be 
noticed. The obtained value is similar to the one characteristic 
for the basic electrolyte. A formation of the composite grains 
(consisting of a single grain and the generated shell) can be 
observed in the sample but the continuous percolation path is 
not yet achieved. A particular filler concentration can be found 
(dependent on the grain size, geometry and the shell thickness) 
for which the conductivity of the system increases abruptly. The 
observed maximum should achieve the values typical for the 
pure amorphous phase. For higher but still moderate amounts of 
the added grains a slow drop of the conductivity value should be 
observed due to the dilution effect. This situation is related to the 
fact that for the composition for which the maximum of conduc- 
tivity is achieved all the polymeric material is already belonging 
to the grain shells. The excessive filler amount cannot, thus, 
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lead to any further improvement of electrical parameters. An 
abrupt depletion is observed at the negative percolation point. 
This phenomenon can be attributed to the decay of continuity of 
the polymer matrix. The charge carrier transport is then blocked 
by the inert insulating filler grains. This phenomenon which is 
predicted by the model can be hardly observed for the experi- 
mental systems while the volume fraction of the filler needed is 
about 0.85. 


1.1. Mesoscale models of conductivity 


First mesoscale model applied for the calculations of conduc- 
tivity in the composite polymeric electrolytes was the Effective 
Medium Theory [7,8]. It was based on the Maxwell—Garnett 
mixing rule for the calculations of effective conductivity of 
the composite unit (grain+ shell) [9] and Landauer [10] and 
Bruggeman [11] works modified by Nan and Smith [12,13] and 
Nakamura [14] applied to the calculation of the final conduc- 
tivity value of the composite for different percolation regimes. 
The main disadvantage of this successful approach is the need 
of an appropriate mixing rule for each studied geometry of the 
composite including such interesting phenomena as shell inho- 
mogeneity, grain size and shape distribution. More versatile but 
on the other hand more complicated approaches to the simulation 
of the composites conductivity are based on the Finite Element 
[15,16] Finite Gradient [17,18] or Random Resistor Network 
(RRN). All three are based on the same concept of building of 
a virtual image of the simulated sample in the form of a two- 
or three-dimensional network of elements representing different 
fragments of the material studied. The main difference lies in the 
used algorithm of the network solving dependent on its varying 
creation rules. 

RRN modeling is a widely recognized approach applied to 
different properties of heterogeneous materials consisting of 
regions of various properties. Apart from the electrical phenom- 
ena, magnetic properties of materials [19] or even the studies 
of sample destruction under critical mechanical stress [20] or 
geophysical studies of earthquakes [21] can be successfully per- 
formed by the means of this computational tool. When applied 
to electrical properties of the heterogeneous systems the method 
relates the value and type of a particular element to local prop- 
erties of the sample. Resistivities, capacitances, parallel R-CPE 
connections can be used being adequate to conductors, insulators 
and other types of materials. Most of already published works 
(e.g. [22—24]) deal with very simple so-called binary systems 
being the image of the two-compound system. Usually, one of 
the components of the mixture is an insulator while the other 
one is a conductor. Electrical parameters of the branches (net- 
work elements) have only two valid values correlated to either 
a path (resistor) or an insulator (break). For such models an 
abrupt percolation threshold can be observed for a given com- 
position of the sample represented by the particular number of 
either conducting (inserting conductors into insulating sample) 
or insulating (breaking the already conducting system) branches 
in the complete network. 

Different more complicated systems were also studied by 
the means of the RRN approach [25-29]. Wang et al. [30] 


predict conductivity in composites consisting of the particles 
of metal and solid electrolyte. The simulated system con- 
sists of units of two different types to which either ionic or 
electronic conductivity can be attributed. To represent differ- 
ent types of materials present in the sample the authors use 
three resistive elements (metal-metal, metal—electrolyte and 
electrolyte—electrolyte connections) and two capacitors (geo- 
metrical capacity of the electrolyte and the double layer capacity 
of the electrolyte—metal junction). The phenomenon of double 
layer creation together with that of respective capacity presence 
is revealed at grain boundaries of different materials. An increase 
of average sample conductivity can be observed for composi- 
tions in which the percolation of metallic grains is observed. 

A model of a system with similar geometry such as composite 
polymeric electrolytes with nonconductive filler is presented in 
[14]. Inorganic composites based on Lil/Al203 and CuCl/Al.03 
are studied. The authors assume that the highly conductive 
interfacial layer is built on the grain—matrix border. The main 
difference between the studied system and the polymeric sys- 
tems lies in the fact that the insulating grains of the inorganic 
composite are not fully randomly located in the volume of the 
real sample. In the presented model randomly located (thus, the 
model is not in agreement with the real geometry) nodes (grains) 
are insulators. Sites not filled with the grain and, thus, belonging 
to the matrix exhibit ionic conductivity. The resistivity value for 
a branch is calculated on the basis of the status of its end nodes. 
If both belong to the grains the branch is broken (infinite resis- 
tivity) and if both belong to the matrix, the branch is the pristine 
conductor. In the case of a branch with ends placed in the nodes 
attributed to two different materials a high conductivity value is 
attributed. Finally, we must conclude that models presented in 
the literature from the point of view of the application in com- 
posite polymeric electrolytes are strongly oversimplified and, 
thus, cannot be applied without further elaboration. The main 
stress must be put on the finite thickness of the shell and its 
heterogeneity. 


1.2. The proposed model 


The proposed model is more sophisticated. It must be stressed 
that realization of so complicated calculations was possible only 
because of the increasing availability of high speed processors 
and large capacity memories in last few years. Finally, it is pos- 
sible to deal with a 300° matrix filled with impedance data for 
the virtual sample represented as complex numbers. The initial 
part of the algorithm is related to the procedure of creation of 
a three-dimensional resistor network being a representation of 
morphology and phase structure of a sample of the simulated 
composition. At the beginning a virtual sample is built of the 
empty cubes. Then the location of spherical grains within the 
matrix is started. The procedure is based on random attributing 
a set of coordinates to the grain center. If a particular simulation 
assumes that grains are not identical a radius of a particular grain 
is also randomly chosen from a given range. Next, all cubic clus- 
ters belonging to the grain are marked using the following rule: 
the cluster belongs to the grain if the center of the cube lies inside 
it. This procedure is repeated until the assumed volume fraction 
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of the filler is achieved. As the ceramic material is not elastic 
and grains cannot overlay, each spatial unit can be attributed to 
only one grain. If a conflict of spatial restrictions is observed, 
during the next grain location an additional procedure changes 
the grain center coordinates and the check is performed once 
again. The algorithm also takes into consideration the insertion 
of grain fragments in the border area of the sample to observe 
the condition that analyzed fragment is a valid representation of 
a larger sample. 

The residual part of the virtual sample (not belonging to 
grains) is then attributed to the shells (respecting the assumed 
shell thickness) and polymer electrolyte matrix (all cells belong- 
ing to neither the grain nor the shell). The values of characteristic 
impedances are attributed to all unitary cells on the basis of the 
assumed material parameters of all phases present in the sys- 
tem. The equivalent circuit used for the representation of each 
cell was a parallel combination of a resistor and a capacitor. 
At the end of this stage the virtual material is described as a 
three-dimensional network of impedances Zyyz. Each cubic unit 
is electrically connected to its six neighbors (only the connec- 
tions by the common walls are taken into consideration). The 
cell-oriented structure is then converted into a branch-oriented 
one with consequent location of the network nodes in the centers 
of the units. Each branch is, thus, characterized by an effective 
impedance value calculated as a serial connection of fractions 
originating from both neighboring cells. 

A branch-oriented three-dimensional network obtained in 
this way can be solved in a various ways including both ana- 
lytical and numerical approaches. For a matrix of the used size 
only numerical iteration procedures can be used. The applied 
algorithm uses a modification of the method introduced by 
Kirkpatrick [24] to find mean parameters of the heterogeneous 
network. To find the impedance of the sample the overall current 
must be calculated for a given test potential difference applied to 
the network. The current value is equal to the one found for any 
of the perpendicular layers of nodes in the sample. This value 
is found as a sum of currents in all branches located between 
nodes belonging to two neighboring layers. This observation is 
valid only in the stationary condition when the first Kirchoff’s 
law is obeyed for each node (1). 


yh =0 (1) 


The network balancing is obtained by the attribution of proper 
potential values to all nodes in the network. Initial potentials val- 
ues must be found as a starting point for the calculations. Al the 
nodes located in the first and the last layer and, thus, attributed to 
the electrodes are fixed to potential values equal 0 and U respec- 
tively. The values of the intra-layer branch impedances in these 
layers are equal zero and the potentials are time-independent. 
Between them monotonic change of the potential is maintained 
by creating a set of equipotential layers of nodes. The value 
of the potential change between two neighboring layers is at 
this step constant. The initial potential distribution obtained in 
this way is later modified by the iteration procedure. In each 
of the iteration steps the value of each node potential is deter- 
mined. A sum of currents for a node is calculated on the basis 


of the actual potential of the node, potentials of the neighboring 
nodes and the impedances of the braches located between them. 
Consequently, the potential of the node is changed in a manner 
allowing the fulfillment of the equation [1]. If it is not possible 
(as the value of the new node potential must lie in the range of 
the neighboring nodes voltages) the modified value is chosen to 
balance the node in the best possible manner. After performing 
the potential shifts for all nodes the whole procedure is repeated 
until the voltage distribution is stabilized. Consequently, the sta- 
tionary state of the branch currents for all of the nodes should 
be reached in this way. The main disadvantage of the procedure 
described above lies in the fact that a quick stabilization of the 
local voltage distribution is reached, unfortunately, in the way 
which prevents global current stabilization. Therefore, an addi- 
tional sub-procedure was created. It takes into consideration the 
interlayer currents forcing an additional modification of the node 
potentials in the beginning of the iteration process. 

The total equalization of the currents in a 100° network is 
obtained in 10-105 iterations. This value is in the upper limit 
range taking into consideration the acceptable calculation time. 
Thus, for the used size of the matrix equal 300° there was a 
need of the change of the calculations termination criteria. The 
analysis of current changes for the following iteration steps in 
connection with an extrapolation procedure to obtain the final 
current value of the balanced network leads to at least one order 
of magnitude decrease in the needed calculation time allowing 
for the calculations of the large samples. We have proved by the 
test simulations that this modification does not affect the final 
conductivity value. In general the used model is a development 
of the dc model [31-33] previously described by us. 


2. Experimental 


OPEN MOSIX working on Slackware GNU/Linux was used 
to manage the computational cluster. C++ was used for code 
generation. The GNU C compiler was used to obtain the exe- 
cutables. The calculation times were in the range of a week for a 
set of data starting from 0 vol.% up to about 35 vol.% of the non- 
conductive additive. Electrical parameters for phases present in 
the system were given as follows: polymer matrix: ø = 1 x 1078, 
k=4; shell: o=1 x 1074, k=8; filler: c= 1 x 107!?, k= 10. For 
some simulations (namely Fig. 4) modified values of the dielec- 
tric constants were applied to achieve a clarity of the phenomena 
described. 


3. Results and discussion 


The first interesting point of the study is the analysis of phase 
composition of the simulated electrolyte. A composite with grain 
diameter equal to 7 um and shell thickness equal to 5 um was 
chosen as a typical representative of the wide range of the mate- 
rials of our interest. The concentrations of the particular phases 
present in the virtual sample are depicted in Fig. 1 as a function 
of the assumed amount of the added ceramic filler. Focusing on 
the shell amount one can observe a maximum located at about 
15 vol.% of the filler. This observation is in good agreement with 
the maximum of the dc conductivity value (see Fig. 2). A second 
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Fig. 1. Phase composition of the simulated samples as a function of the assumed 
filler concentration: (@) grains, (MH) grain shells, (A) polymeric electrolyte 
matrix. 


characteristic point is located at about 25 vol.% of the filler and 
can be described as an increase of the slope of the curve (quicker 
decrease of shell concentration starts here). When comparing 
with the curve related to the pristine polymeric matrix one can 
observe that this threshold is related to the total conversion of the 
matrix into shells (when the matrix concentration reaches zero). 
The conductivity plot also changes here its slope to more nega- 
tive. This set of dependencies together with both linear decrease 
of conductivity and shell volume fraction above this point shows 
clearly that in this filler concentration range (25-30 vol.%) the 
model predicts conductivity changes related to the dilution of 
already fully converted polymeric material by the excess of the 
filler. 
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Fig. 2. Simulated dc conductivity of the set of samples as a function of the filler 
amount. Values of conductivity calculated from the ac model for frequency 
equal 0. 


The comparison of the shell phase amount and virtual sam- 
ple conductivity as a function of the filler contents for samples 
characterized by various grain diameter and constant shell thick- 
ness is shown in the Fig. 3a and b, respectively. One can easily 
observe that the shape of the curve describing the phase composi- 
tion (a) is strongly dependent on the assumed grain diameter. An 
increase of this parameter leads, firstly, to a shift of the maximum 
location to higher filler amounts. Secondly, the maximal value is 
lower for higher grain diameters. Additionally, the slope of the 
increasing branch of the curve is lower and the maximum is more 
flat. For very high grain concentrations all curves overlap (about 
30 vol.%) showing that above that point only a dilution of the 
shells originating from the already converted polymer matrix by 
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Fig. 3. (a) A set of the shell volume fraction dependencies on the filler amount 
for samples with various grain diameter and constant shell thickness (5 um), (b) 
a set of the conductivity dependencies on the filler amount for the same samples: 
(4) d=3 um, (M) d=5 pm, (A) d=7 pm, (@) d=9 pm. 
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Fig. 4. A comparison of the simulated high frequency dielectric response of 
the 10 vol% composite sample: (@) grain diameter d=5 um, shell thickness 
t=3 um, (O) grain diameter d= 5 um, shell thickness t=7 pm. 


the excessive filler grains can be observed. Similar dependencies 
can be observed on the conductivity curve (Fig. 3b) showing that, 
in general, the conductivity value changes follow the changes in 
the amount of the shell phase. 

Fig. 4 shows the simulated high frequency dielectric response 
for two virtual samples of the composite electrolyte containing 
10 vol.% of the filler characterized by the grain diameter equal 
to 5 um. The difference lies in the applied thickness of the amor- 
phous shell formed on the grain surface. In the first case (filled 
dots) it is equal to 3 wm in the latter (empty dots) to 7 wm. 
One can easily observe that the change of thickness leading, in 
consequence, to the change of the network units allocated to 
shells results not only in the change of the impedance value (a 
decrease of about four times for the more modified sample) but 
also to the change of the conductivity mechanism. In the first 
case the amorphous shells do not form the continuous perco- 
lation paths yet. A conductivity value is intermediate between 
the pristine electrolyte and the maximum possible. Both mate- 
rials (pristine polymer and amorphous shells) are characterized 
by different impedance properties leading to a different time 
constant of the respective Nyquist semi-circle (t = RC). Two 
merged semi-circles observed can be attributed to the parts of 
current flowing through shell phase and the pristine polymer. In 
the case of the assumed stronger interaction of the filler with the 
polymer 10 vol.% composition is characterized with almost full 
conversion of the polymeric material. Thus, the only single time 
constant can be observed in the spectra giving as a result only 
one semicircle. 

To compare the calculated ac high frequency dielectric 
response of composite virtual samples containing various 
amount of the filler (and constant filler diameter equal to 
d=7 um and shell thickness equal to t= 5 pm) one must apply 
a modified Nyquist representation. Because of the very broad 
range of the impedance value changes a logarithmic scale was 
used to present the data for the unmodified sample and for the 
composite electrolytes on one plot (Fig. 5). It can be observed 
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Fig. 5. A set of high frequency simulated impedance spectra calculated for a 
virtual sample characterized by d=7 um and t=5 wm. Filler fraction equal to: 
(@) 0 vol.%, (MI) 1 vol.%, (A) 2 vol.%, (@) 5 vol.%, (0) 10 vol.%, (O) 20, ($) 
30 vol.%. 


that most significant difference in the values of the impedance 
can be found for virtual samples with the filler amount between 
0 and 10 vol.%. Later, the changes of both imaginary and real 
part of the impedance are smaller. A monotonically decreasing 
tendency is observed between 0 and 20 vol.% of the filler. For 
30 vol.% sample a slight increase of the impedance is observed as 
the sample is located on the negatively sloped branch of the curve 
depicting the conductivity-composition dependence. Addition- 
ally, one can observe that the shape of the curves is identical 
for the pristine and highly modified electrolyte. In contrast, the 
curves located near the abrupt conductivity change composition 
(2 and 5 vol.%) exhibit a slightly different shape. This phe- 
nomenon can be also attributed to the mechanism of conductivity 
based on simulataneous charge transportation in the electrolyte 
bulk and in highly conductive shells for these samples. 

A more detailed difference can be observed using the linear 
(classical) Nyquist plot. In this case data for both pristine and 
modified systems cannot be clearly depicted on one graph. Thus, 
a comparison of moderately and highly doped samples is pre- 
sented in Fig. 6. The set of impedance spectra was calculated 
for samples of various filler concentration and for grain diam- 
eter d=3 um and shell thickness t=5 pm. The assumed values 
allow us to observe the changes of the curve shape for rela- 
tively low filler concentrations as the increase of the shell phase 
amount is quick (compare with Fig. 3a). Here one can perfectly 
see that while most of the tested compositions lead to one semi- 
circle image the lowest one (1%) gives a strongly deformed 
spectrum in which the contribution of the matrix and shells can 
be found. In this sample the total conductivity is also a com- 
bination of bulk and shell charge transport, both characterized 
by different types of material involved (semi-crystalline matrix 
and amorphous shell). Two partially separated time constants 
can be, thus, attributed to the mechanism of conductivity based 
on simultaneous charge transportation in the electrolyte bulk and 
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Fig. 6. A set of high frequency simulated impedance spectra calculated for a 
virtual sample characterized by d=3 wm and t=5 wm. Filler fraction equal to: 


(¢) 1, (Œ) 2, (0) 3, (O) 5, (@) 10, (I) 20, (4) 30 vol.%. 


in highly conductive shell. Also an interesting phenomenon can 
be observed for 2 and 30 vol.% samples which both give the 
almost identical impedance spectra in the form of a depressed 
semi-circle. This observation is also valid for the 3 and 20 vol.% 
pair and for the 5 and 10 vol.% pair of samples. In the later case, 
the obtained semi-circles are only slightly deformed as the sam- 
ple compositions are located around the conductivity maximum. 
The case of the 3 and 20 vol.% pair is between of the middle of 
the two terminal pairs. Finally, it must be stressed that for moder- 
ate and strongly modified samples only one time constant can be 
observed because only one charge transport mechanism related 
to the amorphous shells is present in the sample. The similarity 
of the pairs of images found below and above the maximum of 
the conductivity is related to the oversimplification of the model 
based on the assumption of the shell homogeneity. As was pre- 
viously proved, for dc calculations [33] this approach produces 
a partially false image of the samples containing large amounts 
of the dispersed filler. 

Finally, the problem of the method used for the calculation 
of the dc conductivity value from the ac impedance spectra can 
be analyzed. A reference method is the extrapolation of the Z’ 
values to the frequency equal zero. This method gives values 
identical with the ones given by the direct application of the 
previously tested dc RRN approach [31-33]. For comparison, 
the values calculated from the diameter of the ac semi-circle 
and from its zero point were plotted together (Fig. 7). For both 
semi-circle based methods the values were calculated not taking 
into consideration the spectrum depression and time constants 
separation (treating data as a one semi-circle). One can see that 
the shapes of the obtained curves are identical and the obtained 
values are very similar. A relative error plot was depicted in 
Fig. 8. One can see that both methods work in the whole range 
of filler amount with the relative error value lower than 5%. The 
method based on the zero point is much more adequate produc- 
ing a few times lower relative error values. Both curves have 
similar shapes with the minimum located near the maximum 
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Fig. 7. Simulated dc conductivity values found from the ac spectra of the virtual 
sample as a function of filler concentration. Calculations performed in various 
manners: (Q) extrapolation of the Z’ to w=0, (O) diameter of the semi-circle, 
(A) zero point of the semi-circle. 


of conductivity. This observation can be related to the fact that 
for these samples the obtained spectra are much closer to the 
ones predicted for the parallel RCPE connection. The discrep- 
ancy of the error values is more pronounced for the high filler 
amount where the obtained spectra form a depressed semi-circle 
showing that in this case the zero-point method produces results 
closer to the expected ones. A smaller difference observed for 
the slightly modified samples shows that both methods are not 
ideal for a merged two semi-circle spectrum characteristic for the 
regime in which mechanism of conductivity is based on simulta- 
neous charge transportation in the electrolyte bulk and in highly 
conductive shells. 
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Fig. 8. Relative errors of the simulated dc conductivity values found from the ac 
spectra of the virtual sample as a function of filler concentration. Calculations 
performed in various manners (extrapolation of the Z’ to œ =0 as a reference): 
(0O) diameter of the semi-circle, (A) zero point of the semi-circle. 
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4. Conclusions 


The Random Resistor Network approach was successfully 
applied to the simulation of ac electrical properties of composite 
polymeric electrolytes. The obtained data are in good agreement 
with the previously published data coming both from impedance 
spectroscopy measurements and the Effective Medium Theory 
simulations. The observed changes of the computed conductivity 
value follow, in general, the observed phase structure modifi- 
cation leading to an increase of the highly conductive phase 
content. For a small amount of grains added one can observe 
an increase of the conductivity value which is also present for 
the virtual systems studied. The slope of the curve and the posi- 
tion of the maximum vary with the changes of grain size and 
the shell thickness. The behavior observed can be attributed to 
the increase of the highly conductive shell contents. Addition- 
ally, a second time constant can be found in some impedance 
spectra of the slightly modified samples for filler amount val- 
ues around the abrupt conductivity change composition. The 
conductivity maximum can be attributed to the total conversion 
of the sample volume (disregarding the units belonging to the 
grains) into highly conductive amorphous shells. No more spa- 
tial units can be attributed to the pristine polymeric matrix. In 
the ac image this range can be described as the well defined 
single semi-arc type spectrum. The following conductivity drop 
is related to a decrease in the number of clusters containing the 
polymeric material due to the dilution effect of the excessive 
filler grains. In the ac spectra this phenomenon can be attributed 
to the depression of the previously observed single semi-arc. 

The proposed model is very versatile and can be easily mod- 
ified to cover systems more complicated and, thus, closer to 
reality. The main possible modification (already realized for the 
dc approach) is related to the introduction of the defined con- 
ductivity distribution within the highly conductive amorphous 
shell. This operation should allow taking into the consideration 
the stiffening phenomenon observed for composites containing 
ceramic fillers. The exchange of nonconductive grains for con- 
ductive ones is also possible and allows a study of so-called 
mixed-phase composites together with their grain boundary 
properties. 
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